Congreso SP(2)
Presentación principales resultados. Parte 2.
Cargar datos y exploración
Cargar paquetes estadísticos y gestión de BBDD
Code
# invisible("Only run from Ubuntu")
# if (!(Sys.getenv("RSTUDIO_SESSION_TYPE") == "server" || file.exists("/.dockerenv"))) {
# if(Sys.info()["sysname"]!="Windows"){
# Sys.setenv(RETICULATE_PYTHON = "/home/fondecytacc/.pyenv/versions/3.11.5/bin/python")
# }
# }
#clean enviroment
rm(list = ls()); gc()
time_before_dedup2<-Sys.time()
if(!require(reticulate)){install.packages("reticulate")}
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
# Busca .mamba_root/envs/py311/python.exe desde getwd() hacia padres
find_python_rel <- function(start = getwd(),
rel = file.path(".mamba_root","envs","py311","python.exe")) {
cur <- normalizePath(start, winslash = "/", mustWork = FALSE)
repeat {
cand <- normalizePath(file.path(cur, rel), winslash = "/", mustWork = FALSE)
if (file.exists(cand)) return(cand)
parent <- dirname(cur)
if (identical(parent, cur)) return(NA_character_) # llegó a la raíz
cur <- parent
}
}
# --- Bootstrap reticulate con ruta relativa a getwd() ---
if(Sys.info()["sysname"]!="Windows"){
#Sys.setenv(RETICULATE_PYTHON = "usr/bin/python3")
reticulate::py_config()
} else {
py <- find_python_rel()
if (is.na(py)) {
stop("No se encontró Python relativo a getwd() (buscando '.mamba_root/envs/py311/python.exe').\n",
"Directorio actual: ", getwd())
}
# Forzar ese intérprete
Sys.unsetenv(c("RETICULATE_CONDAENV","RETICULATE_PYTHON_FALLBACK"))
Sys.setenv(RETICULATE_PYTHON = py)
reticulate::use_python(py, required=T)
reticulate::py_config() # verificación
}
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
resolve_doc_dir <- function() {
# 1) Interactivo en RStudio (mejor para desarrollo)
if (interactive() && requireNamespace("rstudioapi", quietly = TRUE) && rstudioapi::isAvailable()) {
path <- tryCatch(rstudioapi::documentPath(), error = function(e) NULL)
if (!is.null(path) && nzchar(path)) return(normalizePath(dirname(path)))
}
# 2) Render / knitr (Quarto)
in_knitr <- tryCatch(knitr::current_input(), error = function(e) NULL)
if (!is.null(in_knitr) && nzchar(in_knitr)) return(normalizePath(dirname(in_knitr)))
# 3) Si nada funciona, devolver error claro (no fallback silencioso)
stop("No se pudo determinar la carpeta del documento. Asegura que estás renderizando con Quarto o usando RStudio.")
}
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
if(grepl("ubuntu",pak::system_r_platform())){
if (rstudioapi::isAvailable()) {
# Code that needs RStudio
wdpath<- dirname(dirname(rstudioapi::documentPath()))
} else {
# Code for non-RStudio environments (e.g., command line, server)
# This part should use a portable method like knitr::current_input() or getwd()
wdpath <- dirname(resolve_doc_dir())
}
}
load(paste0(wdpath,"/congresosp.RData")) used (Mb) gc trigger (Mb) max used (Mb)
Ncells 602966 32.3 1320207 70.6 835224 44.7
Vcells 1138441 8.7 8388608 64.0 1930478 14.8
python: /home/andres/.cache/R/reticulate/uv/cache/archive-v0/9jodbrlaYfTdFRuQinKW1/bin/python3
libpython: /home/andres/.cache/R/reticulate/uv/python/cpython-3.11.13-linux-x86_64-gnu/lib/libpython3.11.so
pythonhome: /home/andres/.cache/R/reticulate/uv/cache/archive-v0/9jodbrlaYfTdFRuQinKW1:/home/andres/.cache/R/reticulate/uv/cache/archive-v0/9jodbrlaYfTdFRuQinKW1
virtualenv: /home/andres/.cache/R/reticulate/uv/cache/archive-v0/9jodbrlaYfTdFRuQinKW1/bin/activate_this.py
version: 3.11.13 (main, Jul 11 2025, 22:43:55) [Clang 20.1.4 ]
numpy: /home/andres/.cache/R/reticulate/uv/cache/archive-v0/9jodbrlaYfTdFRuQinKW1/lib/python3.11/site-packages/numpy
numpy_version: 2.3.3
NOTE: Python version was forced by py_require()Code
#https://github.com/rstudio/renv/issues/544
#renv falls back to copying rather than symlinking, which is evidently very slow in this configuration.
renv::settings$use.cache(FALSE)
#only use explicit dependencies (in DESCRIPTION)
renv::settings$snapshot.type("implicit")
#check if rstools is installed
if(Sys.info()["sysname"]=="Windows"){
try(installr::install.Rtools(check_r_update=F))
}
check_quarto_version <- function(required = "1.7.29", comparator = c("ge","gt","le","lt","eq")) {
comparator <- match.arg(comparator)
current <- package_version(paste(unlist(quarto::quarto_version()), collapse = "."))
req <- package_version(required)
ok <- switch(comparator,
ge = current >= req,
gt = current > req,
le = current <= req,
lt = current < req,
eq = current == req)
if (!ok) {
stop(sprintf("Quarto version check failed: need %s %s (installed: %s).",
comparator, required, current), call. = FALSE)
}
invisible(TRUE)
}
check_quarto_version("1.7.29", "ge")
#change repository to CL
local({
r <- getOption("repos")
r["CRAN"] <- "https://cran.dcc.uchile.cl/"
options(repos=r)
})
if(!require(pacman)){install.packages("pacman");require(pacman)}Code
if(!require(pak)){install.packages("pak");require(pak)}Code
pacman::p_unlock(lib.loc = .libPaths()) #para no tener problemas reinstalando paquetesCode
if(Sys.info()["sysname"]=="Windows"){
if (getRversion() != "4.4.1") { stop("Requires R version 4.4.1; Actual: ", getRversion()) }
}
#check docker
check_docker_running <- function() {
# Try running 'docker info' to check if Docker is running
system("docker info", intern = TRUE, ignore.stderr = TRUE)
}
if(Sys.info()["sysname"]=="Windows"){
install_docker <- function() {
# Open the Docker Desktop download page in the browser for installation
browseURL("https://www.docker.com/products/docker-desktop")
}
# Main logic
if (inherits(try(check_docker_running(), silent = TRUE), "try-error")) {
liftr::install_docker()
} else {
message("Docker is running.")
}
}
# if(Sys.info()["sysname"]!="Windows"){
# system("sudo apt-get update && sudo apt-get install -y libudunits2-dev libgdal-dev libgeos-dev libproj-dev")
# system("sudo apt-get install -y libudunits2-dev libgdal-dev libgeos-dev libproj-dev")
# }
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
#PACKAGES#######################################################################
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
unlink("*_cache", recursive=T)
# Load or install required packages using pacman::p_load
# This handles installation, loading, and avoids duplicates
# Note: INLA requires separate installation from its dedicated repository.
pacman::p_load(
# Package Management
pacman, # Streamlines loading and installing packages simultaneously
# Core Data & Visualization (Part of tidyverse)
tidyverse, # Collection of packages for data manipulation, visualization, and more
dplyr, # Grammar of data manipulation (e.g., filter, mutate, group_by)
tidyr, # Tools for tidying data, reshaping (e.g., pivot_wider, pivot_longer)
ggplot2, # Grammar of graphics for powerful data plotting
stringr, # Consistent wrappers for common string manipulation tasks
purrr, # Tools for working with functions and vectors/lists (e.g., mapping)
viridisLite, # Provides the perceptually uniform 'viridis' color palettes
scales, # Graphical scales for ggplot2 (e.g., custom breaks, formatting labels for axes)
# Survival & Epidemiology
mexhaz, # Flexible parametric hazard regression models for survival analysis
relsurv, # Core package for **Relative Survival Analysis** in population-based studies
survminer, # Extends 'survival' for high-quality Kaplan-Meier plots and survival tables
popEpi, # Tools for standardized mortality ratios (SMRs) and population-level survival measures
epitools, # Epidemiological tools for data analysis (e.g., rate ratios, exact confidence intervals)
# Spatial Analysis
sf, # Simple Features: for handling and analyzing **vector-based spatial data** (points, lines, polygons)
spdep, # Spatial dependence: for analyzing spatial autocorrelation and spatial regression models
geodata, # Access to global spatial data (e.g., administrative boundaries, climate data)
INLA, # Integrated Nested Laplace Approximations (Bayesian spatial modeling) - install separately
# Data Cleaning & Reporting
janitor, # Simple tools for examining and cleaning dirty data (e.g., cleaning column names)
tableone, # Creates "Table 1" summaries for descriptive statistics (baseline characteristics)
kableExtra, # Enhances 'knitr::kable()' for complex and professional R Markdown tables
pander, # A tool for rendering R objects into Pandoc markdown formats
DT, # Creates interactive HTML data tables (useful for detailed tables in R Markdown/Quarto)
rio, # Simplifies data import/export with a consistent interface ('import()')
quarto, # For rendering the R Markdown/Quarto poster documents (.qmd files)
# Statistical & Computing Tools
biostat3, # Functions and datasets for biostatistics teaching and research
coin, # Conditional inference procedures for robust hypothesis testing (e.g., permutation tests)
metafor, # For heterogeneity testing and general **meta-analysis** methods
parallel, # Base R package for **parallel computing** (useful for speeding up bootstrap and resampling)
cowplot, # Streamlined plot theme and functions for arranging multiple ggplot2 plots
grid, # Base R package for low-level graphics (used indirectly for plot arrangement)
ellmer, # Commonly used for likelihood estimation or generalized mixed-effects models (LMM/GLMM)
# Development & Reproducibility Tools
devtools, # Tools to make package development easier (often used to install from GitHub)
installr, # Functions for self-updating R and R packages
liftr, # A framework for creating and managing reproducible research projects
pak, # A fast, dependency-aware package installer (alternative to install.packages)
install = TRUE
)
# ----------------------------------------------------------------------
# 4. BPMN from GitHub (not on CRAN, so install via devtools if missing)
# ----------------------------------------------------------------------
if (!requireNamespace("bpmn", quietly = TRUE)) {
devtools::install_github("bergant/bpmn")
}
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
#FUNCTIONS######################################################################
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#
# NO MORE DEBUGS
options(error = NULL) # si antes tenías options(error = recover) o browser)
options(browserNLdisabled = FALSE)
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#
#NAs are replaced with "" in knitr kable
options(knitr.kable.NA = '')
pander::panderOptions('big.mark', ',')
pander::panderOptions('decimal.mark', '.')
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#
#to format rows in bold
format_cells <- function(df, rows ,cols, value = c("italics", "bold", "strikethrough")){
# select the correct markup
# one * for italics, two ** for bold
map <- setNames(c("*", "**", "~~"), c("italics", "bold", "strikethrough"))
markup <- map[value]
for (r in rows){
for(c in cols){
# Make sure values are not factors
df[[c]] <- as.character( df[[c]])
# Update formatting
df[r, c] <- ifelse(nchar(df[r, c])==0,"",paste0(markup, gsub(" ", "", df[r, c]), markup))
}
}
return(df)
}
#To produce line breaks in messages and warnings
knitr::knit_hooks$set(
error = function(x, options) {
paste('\n\n<div class="alert alert-danger" style="font-size: small !important;">',
gsub('##', '\n', gsub('^##\ Error', '**Error**', x)),
'</div>', sep = '\n')
},
warning = function(x, options) {
paste('\n\n<div class="alert alert-warning" style="font-size: small !important;">',
gsub('##', '\n', gsub('^##\ Warning:', '**Warning**', x)),
'</div>', sep = '\n')
},
message = function(x, options) {
paste('<div class="message" style="font-size: small !important;">',
gsub('##', '\n', x),
'</div>', sep = '\n')
}
)
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
#CONFIG #######################################################################
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
options(scipen=2) #display numbers rather scientific number
#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_
# Define the function first
#joins these values with semicolons and optionally truncates the result if it exceeds a specified width.
toString2 <- function(x, width = NULL, ...) {
string <- paste(x, collapse = "; ")
if (missing(width) || is.null(width) || width == 0)
return(string)
if (width < 0)
stop("'width' must be positive")
if (nchar(string, type = "w") > width) {
width <- max(6, width)
string <- paste0(substr(string, 1, width - 3), "...")
}
string
}
normalize_txt <- function(x) {
x|>
stringi::stri_trans_general("Latin-ASCII")|>
tolower()|>
trimws()
}
# Custom function for sampling with a seed
sample_n_with_seed <- function(data, size, seed) {
set.seed(seed)
dplyr::sample_n(data, size)
}
# Function to get the most frequent value
most_frequent <- function(x) {
uniq_vals <- unique(x)
freq_vals <- sapply(uniq_vals, function(val) sum(x == val))
most_freq <- uniq_vals[which(freq_vals == max(freq_vals))]
if (length(most_freq) == 1) {
return(most_freq)
} else {
return(NA)
}
}
sum_dates <- function(x){
cbind.data.frame(
min= as.Date(min(unclass(as.Date(x)), na.rm=T), origin = "1970-01-01"),
p001= as.Date(quantile(unclass(as.Date(x)), .001, na.rm=T), origin = "1970-01-01"),
p005= as.Date(quantile(unclass(as.Date(x)), .005, na.rm=T), origin = "1970-01-01"),
p025= as.Date(quantile(unclass(as.Date(x)), .025, na.rm=T), origin = "1970-01-01"),
p25= as.Date(quantile(unclass(as.Date(x)), .25, na.rm=T), origin = "1970-01-01"),
p50= as.Date(quantile(unclass(as.Date(x)), .5, na.rm=T), origin = "1970-01-01"),
p75= as.Date(quantile(unclass(as.Date(x)), .75, na.rm=T), origin = "1970-01-01"),
p975= as.Date(quantile(unclass(as.Date(x)), .975, na.rm=T), origin = "1970-01-01"),
p995= as.Date(quantile(unclass(as.Date(x)), .995, na.rm=T), origin = "1970-01-01"),
p999= as.Date(quantile(unclass(as.Date(x)), .999, na.rm=T), origin = "1970-01-01"),
max= as.Date(max(unclass(as.Date(x)), na.rm=T), origin = "1970-01-01")
)
}
is_stata_ok <- function(x) {
nchar(x) <= 32 & grepl("^[A-Za-z][A-Za-z0-9_]*$", x)
}
to_ascii_lower <- function(x) stringi::stri_trans_general(x, "Latin-ASCII")
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
# ── Tidy RMS ────────────────────────────────────────────────────────
tidy_cph <- function(model) {
if (!inherits(model, "cph")) stop("Model must be an rms::cph object")
# run summary
s <- summary(model)
df_s <- as.data.frame(unclass(s))
df_s$rown <- rownames(s)
# find rows that are "Hazard Ratio"
hr_rows <- trimws(df_s$rown) == "Hazard Ratio"
# build tidy tibble
tibble::tibble(
term = trimws(rownames(s)[which(hr_rows) - 1L]), # previous row = variable name
estimate = df_s$Effect[hr_rows], # Hazard Ratio
conf.low = df_s$`Lower 0.95`[hr_rows],
conf.high= df_s$`Upper 0.95`[hr_rows],
p.value = df_s$P[hr_rows]
)
}
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
# ── SMR Phi ────────────────────────────────────────────────────────
sir_ci_phi_improved <- function(sir_obj, phi, conf.level = 0.95) {
#Método log-normal, the best, dont overestimate, or subestimate variance
# extract totals
total_obs <- sir_obj$observed
total_exp <- sir_obj$expected
# Calculate SEs
theta <- total_obs / total_exp
# Normal approximation, n>20
# Corrected SEs (McCullagh & Nelder, 1989)
# “For ratios of Poisson means (such as SIR or CMR), the appropriate approach is to use multinomial or binomial models conditioned on the total observed.”
# Breslow NE, Day NE. Statistical Methods in Cancer Research, Vol. II (IARC, 1987), §2.2. – Derives the same SE formula and recommends inflating by φ in the presence of heterogeneity.
z <- qnorm(1 - (1 - conf.level)/2)
se_log <- sqrt(phi / total_obs) # Valid formula
# ICs
lci <- theta * exp(-z * se_log)
uci <- theta * exp(z * se_log)
data.frame(
SIR = theta,
CI_low = lci,
CI_high = uci,
phi_used = phi
)
}
smr_print <- function(sir_out, name = NULL) {
# Validate input
if (!is.list(sir_out)) stop("sir_out must be a list or list-like object.")
required <- c("observed", "pyrs", "expected", "sir")
missing_req <- setdiff(required, names(sir_out))
if (length(missing_req)) stop("sir_out is missing required elements: ", paste(missing_req, collapse = ", "))
# Pull components (allow vectors)
obs <- as.numeric(sir_out$observed)
pyrs <- as.numeric(sir_out$pyrs)
expc <- as.numeric(sir_out$expected)
sir <- as.numeric(sir_out$sir)
sir_lo <- if ("sir.lo" %in% names(sir_out)) as.numeric(sir_out$sir.lo) else rep(NA_real_, length(sir))
sir_hi <- if ("sir.hi" %in% names(sir_out)) as.numeric(sir_out$sir.hi) else rep(NA_real_, length(sir))
ear <- if ("EAR" %in% names(sir_out)) as.numeric(sir_out$EAR) else rep(NA_real_, length(sir))
n <- max(length(obs), length(pyrs), length(expc), length(sir), length(sir_lo), length(sir_hi), length(ear))
rep_len <- function(x) if (length(x) == n) x else rep(x, length.out = n)
obs <- rep_len(obs)
pyrs <- rep_len(pyrs)
expc <- rep_len(expc)
sir <- rep_len(sir)
sir_lo <- rep_len(sir_lo)
sir_hi <- rep_len(sir_hi)
ear <- rep_len(ear)
name <- if (is.null(name)) rep(NA_character_, n) else rep_len(as.character(name))
# Format SMR: "x.xx (y.yy–z.zz)" or "x.xx" if bounds missing
smr_fmt <- vapply(seq_len(n), function(i) {
if (!is.na(sir_lo[i]) && !is.na(sir_hi[i])) {
sprintf("%.2f (%.2f–%.2f)", sir[i], sir_lo[i], sir_hi[i])
} else {
sprintf("%.2f", sir[i])
}
}, FUN.VALUE = character(1), USE.NAMES = FALSE)
# Format EAR
ear_fmt <- ifelse(is.na(ear), NA_character_, sprintf("%.2f", ear))
# Build result
out <- data.frame(
total = name,
observed = round(obs, 0),
pyrs = round(pyrs, 0),
expected = round(expc, 0),
SMR = smr_fmt,
EAR = ear_fmt,
stringsAsFactors = FALSE,
row.names = NULL
)
out
}
# ── Plots ────────────────────────────────────────────────────────
theme_sjPlot_manual <- function() {
theme(
plot.title = element_text(size = 16, face = "bold"),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12),
axis.text.y = element_text(size = 12),
panel.background = element_rect(fill = "white"),
panel.grid.major = element_line(color = "gray80"),
panel.grid.minor = element_line(color = "gray90")
)
}
# --- Paleta base continua (tu cont_palette) ---
andal_pal <- function(reverse = FALSE) {
cols <- c(
"#052709", # dark teal
"#5FA9B3", # teal accent
"#A17C6C", # warm mid
"#D2A899", # warm light
"#906F00", # earth gold
"#527700" # olive
)
if (reverse) cols <- rev(cols)
grDevices::colorRampPalette(cols, interpolate = "spline")
}
# Continua (equivalente a scale_fill_viridis_c)
scale_fill_andal_c <- function(..., reverse = FALSE, na.value = "grey85",
guide = "colourbar") {
pal <- andal_pal(reverse)
ggplot2::scale_fill_gradientn(colours = pal(256), na.value = na.value,
guide = guide, ...)
}
scale_color_andal_c <- function(..., reverse = FALSE, na.value = "grey85",
guide = "colourbar") {
pal <- andal_pal(reverse)
ggplot2::scale_colour_gradientn(colours = pal(256), na.value = na.value,
guide = guide, ...)
}
# Binned / escalonada (tipo scale_fill_viridis_b)
scale_fill_andal_b <- function(..., reverse = FALSE, na.value = "grey85",
n.breaks = 6) {
cols <- c(
"#052709", "#5FA9B3", "#A17C6C",
"#D2A899", "#906F00", "#527700"
)
if (reverse) cols <- rev(cols)
ggplot2::scale_fill_stepsn(colors = cols, n.breaks = n.breaks,
na.value = na.value, ...)
}
scale_color_andal_b <- function(..., reverse = FALSE, na.value = "grey85",
n.breaks = 6) {
cols <- c(
"#052709", "#5FA9B3", "#A17C6C",
"#D2A899", "#906F00", "#527700"
)
if (reverse) cols <- rev(cols)
ggplot2::scale_colour_stepsn(colors = cols, n.breaks = n.breaks,
na.value = na.value, ...)
}
# Discreta (por si necesitas categorías)
scale_fill_andal_d <- function(..., reverse = FALSE) {
cols <- c(
"#052709", "#5FA9B3", "#A17C6C",
"#D2A899", "#906F00", "#527700"
)
if (reverse) cols <- rev(cols)
ggplot2::scale_fill_manual(values = cols, ...)
}
scale_color_andal_d <- function(..., reverse = FALSE) {
cols <- c(
"#052709", "#5FA9B3", "#A17C6C",
"#D2A899", "#906F00", "#527700"
)
if (reverse) cols <- rev(cols)
ggplot2::scale_colour_manual(values = cols, ...)
}Estructurar los datos
Code
# Function to clean and prepare region names for matching
# Load the fuzzy matching library
library(stringdist)Code
ref_regions <- tibble::tribble(
~roman_code, ~canon_name, ~aliases,
"XV", "arica y parinacota", c("arica y parinacota","arica","parinacota"),
"I", "tarapaca", c("tarapaca"),
"II", "antofagasta", c("antofagasta"),
"III", "atacama", c("atacama"),
"IV", "coquimbo", c("coquimbo"),
"V", "valparaiso", c("valparaíso","valparaiso"),
"RM", "metropolitana", c("rm","metropolitana","metropolitana de santiago","santiago"),
"VI", "ohiggins", c("o'higgins","ohiggins","libertador general bernardo ohiggins","bernardo ohiggins","libertador ohiggins"),
"VII", "maule", c("maule"),
"VIII","biobio", c("biobio","bio bio","bio-bio","bío-bío"),
"IX", "araucania", c("la araucania","araucania","araucanía"),
"X", "los lagos", c("los lagos"),
"XI", "aysen", c("aysen","aysén"),
"XII", "magallanes y la antartica chilena", c("magallanes","antartica","antártica","magallanes y antartica","magallanes y la antartica chilena"),
"XIV", "los rios", c("los rios","los ríos")
)
# --- Limpieza homogénea (preserva 'los/las') ---
clean_string <- function(s){
s <- tolower(s)
s <- iconv(s, from = "UTF-8", to = "ASCII//TRANSLIT") # quita tildes
s <- gsub("[[:punct:]]", " ", s) # deja espacios
s <- gsub("[[:digit:]]", "", s)
# quitamos conectores pero NO 'los/las' (desambiguación)
s <- gsub("\\b(de|del|y|la|el|region|región)\\b", " ", s)
s <- gsub("\\s+", " ", s)
trimws(s)
}
# expandimos alias y limpiaremos todos con la misma función
ref_alias <- ref_regions |>
tidyr::unnest_longer(aliases) |>
dplyr::mutate(clean_alias = clean_string(aliases),
clean_canon = clean_string(canon_name)) |>
dplyr::distinct(roman_code, canon_name, clean_alias, clean_canon)
# --- Reglas exactas/regex muy específicas (antes de fuzzy) ---
rule_match_code <- function(x_clean){
# pares conflictivos: exigir frase completa
if(grepl("\\blos\\s+lagos\\b", x_clean)) return("X")
if(grepl("\\blos\\s+rios\\b", x_clean)) return("XIV")
# otros casos robustos
if(grepl("\\bbiobio\\b|\\bbio\\s*bio\\b", x_clean)) return("VIII")
if(grepl("\\bo\\s*higgins\\b|\\blibertador\\b|\\bbernardo\\b", x_clean)) return("VI")
if(grepl("\\bmetropolitana\\b|\\brm\\b|\\bsantiago\\b", x_clean)) return("RM")
if(grepl("\\bmagallanes\\b|\\bantartic", x_clean)) return("XII")
if(grepl("\\baraucan", x_clean)) return("IX")
if(grepl("\\bvalparais", x_clean)) return("V")
if(grepl("\\btarapac", x_clean)) return("I")
if(grepl("\\barica\\b|\\bparinacota\\b", x_clean)) return("XV")
if(grepl("\\bantofagasta\\b", x_clean)) return("II")
if(grepl("\\batacama\\b", x_clean)) return("III")
if(grepl("\\bcoquimbo\\b", x_clean)) return("IV")
if(grepl("\\bmaule\\b", x_clean)) return("VII")
if(grepl("\\baysen\\b", x_clean)) return("XI")
return(NA_character_)
}
# --- Matcher seguro ---
match_regions_safe <- function(messy_names, ref_alias,
jw_thresh = 0.12, tie_margin = 0.02){
sapply(messy_names, function(x){
x_clean <- clean_string(x)
# 1) Reglas
code <- rule_match_code(x_clean)
if(!is.na(code)) return(code)
# 2) Igualdad exacta con alias limpios o canónicos
eq_hit <- ref_alias |>
dplyr::filter(x_clean == clean_alias | x_clean == clean_canon)
if(nrow(eq_hit) == 1) return(eq_hit$roman_code[[1]])
if(nrow(eq_hit) > 1) return(NA_character_) # ambigüedad improbable: revisar
# 3) Fuzzy con umbral y desempate
cand <- unique(c(ref_alias$clean_alias, ref_alias$clean_canon))
d <- stringdist::stringdist(x_clean, cand, method = "jw")
i_min <- which.min(d); d_min <- d[i_min]
# si la mejor distancia es mala, no arriesgar
if(is.infinite(d_min) || is.na(d_min) || d_min > jw_thresh) return(NA_character_)
# empate cercano: si la diferencia entre 1º y 2º es pequeña, exigir frase completa o marcar NA
ord <- order(d); if(length(ord) >= 2){
if((d[ord[2]] - d[ord[1]]) < tie_margin) return(NA_character_)
}
cand_str <- cand[i_min]
hit <- ref_alias |>
dplyr::filter(clean_alias == cand_str | clean_canon == cand_str) |>
dplyr::slice(1)
hit$roman_code[[1]]
})
}
# This is your sample data from the query
matched <- match_regions_safe(CONS_C3_mod_c_a2$region_del_centro, ref_alias)
# Add the new codes back to your original data frame
CONS_C3_mod_c_a2$CC1 <- matched
table(CONS_C3_mod_c_a2$region_del_centro, CONS_C3_mod_c_a2$CC1) |>
knitr::kable("markdown", caption="Tabla de contingencia entre nombres originales y códigos asignados")
# Now you can proceed with your lexpand call, including 'reg_res' in the aggre list
SISTRAT_c3_reg <- popEpi::lexpand(
CONS_C3_mod_c_a2,
status = status,
birth = birth_date_rec_joined,
exit = def_date_na,
entry = fecha_ingreso_a_tratamiento,
# 'reg_res' has been REMOVED from the breaks list
breaks = list(
per = seq(2011, 2021, by = 1), #En popEpi::lexpand, los breaks son cortes y generan intervalos [a, b) (cerrado a la izquierda, abierto a la derecha). Cambio en 13-octubre 2025 2025-10-13
age = c(18, 30, 45, 60, Inf)
),
# 'reg_res' has been ADDED to the aggre list
aggre = list(
agegroup = age,
sex = sex,
CC1 = CC1
)
)
# --- 1) Prepare SISTRAT cohort by CC_1 × sex × age ---
# from0to1 is deaths; pyrs are person-years
SISTRAT_agg <- SISTRAT_c3_reg %>%
mutate(
sex_rec = recode(tolower(sex), male = "hombre", female = "mujer"),
adm_age_cat = case_when(
agegroup == 18 ~ "18-29",
agegroup == 30 ~ "30-44",
agegroup == 45 ~ "45-59",
agegroup == 60 ~ "60-78",
TRUE ~ NA_character_
)
) %>%
filter(!is.na(adm_age_cat), !is.na(CC1)) %>%
group_by(CC1, sex_rec, adm_age_cat) %>%
summarise(
Y = sum(from0to1, na.rm = TRUE),
PY = sum(pyrs, na.rm = TRUE),
.groups = "drop"
)
# Cohort totals by CC_1 (for modeling)
SISTRAT_by_CC1 <- SISTRAT_agg %>%
group_by(CC1) %>%
summarise(Y = sum(Y), PY = sum(PY), .groups = "drop")| I | II | III | IV | IX | RM | V | VI | VII | VIII | X | XII | XIV | XV | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| de antofagasta | 0 | 119 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
| de arica y parinacota | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 102 |
| de atacama | 0 | 0 | 22 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
| de coquimbo | 0 | 0 | 0 | 68 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
| de la araucania | 0 | 0 | 0 | 0 | 42 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
| de los lagos | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 81 | 0 | 0 | 0 |
| de los rios | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 56 | 0 |
| de magallanes y la antartica chilena | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 78 | 0 | 0 |
| de tarapaca | 29 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
| de valparaiso | 0 | 0 | 0 | 0 | 0 | 0 | 196 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
| del bio-bio | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 154 | 0 | 0 | 0 | 0 |
| del libertador general bernardo ohiggins | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 82 | 0 | 0 | 0 | 0 | 0 | 0 |
| del maule | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 41 | 0 | 0 | 0 | 0 | 0 |
| metropolitana | 0 | 0 | 0 | 0 | 0 | 321 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
Definir la población de referencia y su mortalidad
Code
# Load the reference mortality rates
proy_ine_reg_group_25_oct_reg_al2015_dput<-
structure(list(reg_res = c("01", "01", "01", "01", "01", "01",
"01", "01", "02", "02", "02", "02", "02", "02", "02", "02", "03",
"03", "03", "03", "03", "03", "03", "03", "04", "04", "04", "04",
"04", "04", "04", "04", "05", "05", "05", "05", "05", "05", "05",
"05", "06", "06", "06", "06", "06", "06", "06", "06", "07", "07",
"07", "07", "07", "07", "07", "07", "08", "08", "08", "08", "08",
"08", "08", "08", "09", "09", "09", "09", "09", "09", "09", "09",
"10", "10", "10", "10", "10", "10", "10", "10", "11", "11", "11",
"11", "11", "11", "11", "11", "12", "12", "12", "12", "12", "12",
"12", "12", "13", "13", "13", "13", "13", "13", "13", "13", "14",
"14", "14", "14", "14", "14", "14", "14", "15", "15", "15", "15",
"15", "15", "15", "15", "16", "16", "16", "16", "16", "16", "16",
"16"), year = c(2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015,
2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015),
sex = c("Male", "Male", "Male", "Male", "Female", "Female",
"Female", "Female", "Male", "Male", "Male", "Male", "Female",
"Female", "Female", "Female", "Male", "Male", "Male", "Male",
"Female", "Female", "Female", "Female", "Male", "Male", "Male",
"Male", "Female", "Female", "Female", "Female", "Male", "Male",
"Male", "Male", "Female", "Female", "Female", "Female", "Male",
"Male", "Male", "Male", "Female", "Female", "Female", "Female",
"Male", "Male", "Male", "Male", "Female", "Female", "Female",
"Female", "Male", "Male", "Male", "Male", "Female", "Female",
"Female", "Female", "Male", "Male", "Male", "Male", "Female",
"Female", "Female", "Female", "Male", "Male", "Male", "Male",
"Female", "Female", "Female", "Female", "Male", "Male", "Male",
"Male", "Female", "Female", "Female", "Female", "Male", "Male",
"Male", "Male", "Female", "Female", "Female", "Female", "Male",
"Male", "Male", "Male", "Female", "Female", "Female", "Female",
"Male", "Male", "Male", "Male", "Female", "Female", "Female",
"Female", "Male", "Male", "Male", "Male", "Female", "Female",
"Female", "Female", "Male", "Male", "Male", "Male", "Female",
"Female", "Female", "Female"), agegroup = c(18, 30, 45, 60,
18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45,
60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30,
45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18,
30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60,
18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45,
60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30,
45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18,
30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60, 18, 30, 45, 60,
18, 30, 45, 60), poblacion = c(34002L, 38653L, 28302L, 15091L,
32200L, 36648L, 27375L, 16273L, 64007L, 72673L, 53668L, 27557L,
60944L, 67667L, 51267L, 30511L, 28880L, 32559L, 28056L, 16330L,
27350L, 30400L, 27570L, 17325L, 71827L, 78583L, 69004L, 45243L,
71480L, 81367L, 72797L, 51981L, 180134L, 186088L, 167910L,
118195L, 172861L, 186938L, 182401L, 142333L, 79763L, 101327L,
93062L, 58682L, 76880L, 99496L, 92575L, 63692L, 93197L, 110899L,
103897L, 68655L, 96053L, 113384L, 106240L, 75650L, 155816L,
164863L, 152554L, 94068L, 155474L, 171681L, 161491L, 111956L,
88721L, 99564L, 90612L, 59966L, 91078L, 103766L, 92400L,
68850L, 77395L, 96278L, 82763L, 48381L, 75995L, 95497L, 80050L,
54458L, 9138L, 12044L, 10142L, 5664L, 8611L, 11949L, 9318L,
5554L, 15950L, 19178L, 16790L, 10861L, 14332L, 18242L, 16390L,
10938L, 760210L, 819336L, 644988L, 385795L, 727290L, 802028L,
699815L, 487007L, 37117L, 38145L, 39261L, 24363L, 37110L,
39544L, 38848L, 28148L, 23852L, 24685L, 19338L, 12961L, 21661L,
24452L, 20876L, 15206L, 41936L, 47039L, 49239L, 33580L, 43939L,
50966L, 51433L, 37434L)), class = c("tbl_df", "tbl", "data.frame"
), row.names = c(NA, -128L))
#mort_2015_reg
mort_2015_reg_dput<-
structure(list(reg_res = c("01", "01", "01", "01", "01", "01",
"01", "01", "01", "01", "02", "02", "02", "02", "02", "02", "02",
"02", "02", "02", "03", "03", "03", "03", "03", "03", "03", "03",
"03", "03", "04", "04", "04", "04", "04", "04", "04", "04", "04",
"04", "05", "05", "05", "05", "05", "05", "05", "05", "05", "05",
"06", "06", "06", "06", "06", "06", "06", "06", "06", "06", "07",
"07", "07", "07", "07", "07", "07", "07", "07", "07", "07", "08",
"08", "08", "08", "08", "08", "08", "08", "08", "08", "08", "09",
"09", "09", "09", "09", "09", "09", "09", "09", "09", "09", "10",
"10", "10", "10", "10", "10", "10", "10", "10", "10", "10", "11",
"11", "11", "11", "11", "11", "11", "11", "11", "11", "12", "12",
"12", "12", "12", "12", "12", "12", "12", "12", "13", "13", "13",
"13", "13", "13", "13", "13", "13", "13", "13", "14", "14", "14",
"14", "14", "14", "14", "14", "14", "14", "14", "15", "15", "15",
"15", "15", "15", "15", "15", "15", "15", "15"), sexo = c(1,
1, 1, 1, 1, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 1, 1,
1, 1, 1, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 1, 1, 1,
1, 1, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 1, 1, 1, 1,
1, 2, 2, 2, 2, 2, 9, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 9, 1, 1, 1,
1, 1, 2, 2, 2, 2, 2, 9, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 9, 1, 1,
1, 1, 1, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 1, 1, 1,
1, 1, 2, 2, 2, 2, 2, 9, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 9, 1, 1,
1, 1, 1, 2, 2, 2, 2, 2, 9), adm_age_cat = c("18-29", "30-44",
"45-59", "60-78", NA, "18-29", "30-44", "45-59", "60-78", NA,
"18-29", "30-44", "45-59", "60-78", NA, "18-29", "30-44", "45-59",
"60-78", NA, "18-29", "30-44", "45-59", "60-78", NA, "18-29",
"30-44", "45-59", "60-78", NA, "18-29", "30-44", "45-59", "60-78",
NA, "18-29", "30-44", "45-59", "60-78", NA, "18-29", "30-44",
"45-59", "60-78", NA, "18-29", "30-44", "45-59", "60-78", NA,
"18-29", "30-44", "45-59", "60-78", NA, "18-29", "30-44", "45-59",
"60-78", NA, "18-29", "30-44", "45-59", "60-78", NA, "18-29",
"30-44", "45-59", "60-78", NA, NA, "18-29", "30-44", "45-59",
"60-78", NA, "18-29", "30-44", "45-59", "60-78", NA, NA, "18-29",
"30-44", "45-59", "60-78", NA, "18-29", "30-44", "45-59", "60-78",
NA, NA, "18-29", "30-44", "45-59", "60-78", NA, "18-29", "30-44",
"45-59", "60-78", NA, NA, "18-29", "30-44", "45-59", "60-78",
NA, "18-29", "30-44", "45-59", "60-78", NA, "18-29", "30-44",
"45-59", "60-78", NA, "18-29", "30-44", "45-59", "60-78", NA,
"18-29", "30-44", "45-59", "60-78", NA, "18-29", "30-44", "45-59",
"60-78", NA, NA, "18-29", "30-44", "45-59", "60-78", NA, "18-29",
"30-44", "45-59", "60-78", NA, NA, "18-29", "30-44", "45-59",
"60-78", NA, "18-29", "30-44", "45-59", "60-78", NA, NA), n = c(40L,
73L, 134L, 316L, 228L, 11L, 31L, 57L, 196L, 325L, 71L, 99L, 338L,
729L, 483L, 28L, 47L, 144L, 477L, 617L, 23L, 69L, 137L, 321L,
326L, 17L, 20L, 73L, 202L, 366L, 83L, 98L, 355L, 895L, 878L,
21L, 61L, 203L, 561L, 1001L, 179L, 271L, 811L, 2531L, 2516L,
45L, 148L, 536L, 1824L, 3413L, 97L, 172L, 449L, 1167L, 1072L,
33L, 95L, 264L, 795L, 1287L, 120L, 184L, 494L, 1475L, 1326L,
43L, 92L, 292L, 978L, 1485L, 2L, 188L, 384L, 1127L, 2850L, 2462L,
71L, 153L, 667L, 1953L, 3057L, 1L, 98L, 206L, 506L, 1359L, 1331L,
33L, 93L, 303L, 938L, 1631L, 4L, 101L, 202L, 507L, 1026L, 952L,
31L, 82L, 229L, 748L, 1228L, 2L, 15L, 22L, 64L, 119L, 99L, 1L,
7L, 21L, 61L, 89L, 18L, 28L, 89L, 247L, 176L, 6L, 16L, 45L, 164L,
244L, 700L, 1169L, 3089L, 7822L, 6991L, 238L, 565L, 1873L, 5714L,
10914L, 8L, 42L, 96L, 215L, 583L, 549L, 14L, 35L, 106L, 383L,
566L, 1L, 33L, 42L, 99L, 290L, 233L, 8L, 23L, 52L, 182L, 282L,
1L), sex_rec = c("hombre", "hombre", "hombre", "hombre", "hombre",
"mujer", "mujer", "mujer", "mujer", "mujer", "hombre", "hombre",
"hombre", "hombre", "hombre", "mujer", "mujer", "mujer", "mujer",
"mujer", "hombre", "hombre", "hombre", "hombre", "hombre", "mujer",
"mujer", "mujer", "mujer", "mujer", "hombre", "hombre", "hombre",
"hombre", "hombre", "mujer", "mujer", "mujer", "mujer", "mujer",
"hombre", "hombre", "hombre", "hombre", "hombre", "mujer", "mujer",
"mujer", "mujer", "mujer", "hombre", "hombre", "hombre", "hombre",
"hombre", "mujer", "mujer", "mujer", "mujer", "mujer", "hombre",
"hombre", "hombre", "hombre", "hombre", "mujer", "mujer", "mujer",
"mujer", "mujer", "mujer", "hombre", "hombre", "hombre", "hombre",
"hombre", "mujer", "mujer", "mujer", "mujer", "mujer", "mujer",
"hombre", "hombre", "hombre", "hombre", "hombre", "mujer", "mujer",
"mujer", "mujer", "mujer", "mujer", "hombre", "hombre", "hombre",
"hombre", "hombre", "mujer", "mujer", "mujer", "mujer", "mujer",
"mujer", "hombre", "hombre", "hombre", "hombre", "hombre", "mujer",
"mujer", "mujer", "mujer", "mujer", "hombre", "hombre", "hombre",
"hombre", "hombre", "mujer", "mujer", "mujer", "mujer", "mujer",
"hombre", "hombre", "hombre", "hombre", "hombre", "mujer", "mujer",
"mujer", "mujer", "mujer", "mujer", "hombre", "hombre", "hombre",
"hombre", "hombre", "mujer", "mujer", "mujer", "mujer", "mujer",
"mujer", "hombre", "hombre", "hombre", "hombre", "hombre", "mujer",
"mujer", "mujer", "mujer", "mujer", "mujer")), row.names = c(NA,
-157L), class = c("tbl_df", "tbl", "data.frame"))
mort_2015_reg_dput<-
dplyr::filter(data.frame(mort_2015_reg_dput),!is.na(adm_age_cat))
# Mapeo de número romano (CC_1) -> código CUT (reg_res)
roman_to_cut <- tibble::tibble(
CC_1 = c("I","II","III","IV","V","VI","VII","VIII","IX","X","XI","XII","RM","XIV","XV","XVI"),
reg_res = c("01","02","03","04","05","06","07","08","09","10","11","12","13","14","15","16")
)
# --- 2) Build background reference rates by CC_1 × sex × age (2015) ---
# Population
pop <- proy_ine_reg_group_25_oct_reg_al2015_dput|>
filter(year == 2015, sex %in% c("Male","Female"))|>
dplyr::mutate(
sex_rec = recode(sex, Male = "hombre", Female = "mujer"),
adm_age_cat = case_when(
agegroup == 18 ~ "18-29",
agegroup == 30 ~ "30-44",
agegroup == 45 ~ "45-59",
agegroup == 60 ~ "60-78",
TRUE ~ NA_character_
),
reg_res = stringr::str_pad(reg_res, width = 2, pad = "0")
)|>
dplyr::filter(!is.na(adm_age_cat))|>
dplyr::select(reg_res, sex_rec, adm_age_cat, poblacion)|>
dplyr::left_join(roman_to_cut, by = "reg_res")|>
dplyr::rename("CC1" = "CC_1")
# Deaths
deaths <- mort_2015_reg_dput|>
dplyr::filter(!is.na(adm_age_cat), sexo %in% c(1,2))|>
dplyr::mutate(
sex_rec = if_else(sexo == 1, "hombre", "mujer"),
reg_res = stringr::str_pad(reg_res, width = 2, pad = "0"))|>
dplyr::select(reg_res, sex_rec, adm_age_cat, n)|>
dplyr::left_join(roman_to_cut, by = "reg_res")|>
dplyr::rename("CC1" = "CC_1")
# Make complete shells to avoid missing strata
shell_cc1 <- expand.grid(
CC1 = unique(roman_to_cut$CC_1),
sex_rec = c("hombre","mujer"),
adm_age_cat= c("18-29","30-44","45-59","60-78"),
stringsAsFactors = FALSE
) |> tibble::as_tibble()
pop_ref <- shell_cc1|>
left_join(pop, by = c("CC1","sex_rec","adm_age_cat"))|>
mutate(poblacion = replace_na(poblacion, 0L))
deaths_ref <- shell_cc1|>
left_join(deaths, by = c("CC1","sex_rec","adm_age_cat"))|>
mutate(n = replace_na(n, 0L))
# Reference rates r_{CC_1,a,s} = D / P
P_df <- pop_ref |> group_by(CC1, sex_rec, adm_age_cat)|> summarise(P = sum(poblacion), .groups = "drop")
D_df <- deaths_ref%>% group_by(CC1, sex_rec, adm_age_cat)|> summarise(D = sum(n), .groups = "drop")
rates_ref <- full_join(D_df, P_df, by = c("CC1","sex_rec","adm_age_cat"))|>
mutate(
D = coalesce(D, 0L),
P = coalesce(P, 0L),
rate = if_else(P > 0, D / P, 0)
)
# --- 3) Expected deaths for the cohort by CC_1: E = Σ PY × rate_ref ---
E_by_CC1 <- SISTRAT_agg|>
left_join(rates_ref, by = c("CC1","sex_rec","adm_age_cat"))|>
mutate(rate = replace_na(rate, 0))|>
group_by(CC1)|>
summarise(
E = sum(PY * rate),
PY = sum(PY),
Y = sum(Y),
.groups = "drop"
)|>
mutate(E = if_else(E == 0, 1e-8, E)) # guard against exact zerosFormatear el mapa y matriz de adyacencia
Armonizamos el mapa con los datos y construimos la matriz de adyacencia para el modelo BYM, fijando a priori conservadores: asumimos una probabilidad baja (≈1%) de que la desviación estándar del efecto aleatorio supere 1 (con esto, asumimos poca variabilidad interregional) no y predefinimos un equilibrio entre el componente espacial estructurado (por vecindad) y el no estructurado (ruido, o particularidad de la región).
Code
# --- 4) Geometry (join by CC_1) and adjacency for BYM ---
chile_sf <- geodata::gadm(country = "CHL", level = 1, path = tempdir()) |>
sf::st_as_sf() |>
dplyr::mutate(CC_1 = as.character(CC_1))|>
# keep only regions we can map (GADM may or may not have XVI; most recent does)
(\(df) {
if (filter(df, !CC_1 %in% roman_to_cut$CC_1) |> nrow()>0) {
message("Note: data does not include X10=Región de Los Lagos; XI11=Región de Aysén del General Carlos Ibáñez del Campo; XVI16=Región de Ñuble")
print(filter(df, !CC_1 %in% roman_to_cut$CC_1) |>pull(CC1))
}
df
})()|>
dplyr::filter(CC_1 %in% roman_to_cut$CC_1)|>
dplyr::arrange(CC_1)|>
dplyr::mutate(region_id = dplyr::row_number())
# Align modeling data to shapes (right_join keeps all regions in the map)
bym_df <- E_by_CC1|>
dplyr::right_join(sf::st_drop_geometry(chile_sf)|> dplyr::select(CC_1, region_id), by = c("CC1"="CC_1"))|>
dplyr::arrange(region_id)|>
dplyr::mutate(
Y = replace_na(Y, 0),
PY = replace_na(PY, 0),
E = replace_na(E, 1e-8) # tiny to allow model to run if truly no expected deaths
)
# Adjacency
nb <- spdep::poly2nb(chile_sf, queen = TRUE)
adj_file <- file.path(tempdir(), "chile_cc1_adj.graph")
spdep::nb2INLA(file = adj_file, nb)
# Priors
hyper <- list(
prec = list(prior = "pc.prec", param = c(1, 0.01)),
phi = list(prior = "pc", param = c(0.5, 2/3))
)Modelo Bayesiano BYM2 por tasas suavizadas, riesgo relativo y muertes en exceso
Obtenemos las tasas por 100 años-persona, riesgo relativo y muertes en exceso suavizadas por vecindad espacial, además de la probabilidad de que el riesgo relativo sea mayor que 1.
Code
# --- 5A) BYM for smoothed death rates (cohort) ---
# Y ~ Poisson(PY × rate_i); fitted values are on scale of rate (per person-year)
fit_rate <- INLA::inla(
Y ~ 1 + f(region_id, model = "bym2", graph = adj_file, scale.model = TRUE, hyper = hyper),
data = bym_df,
family = "poisson",
E = PY, # exposure = person-years
control.predictor = list(compute = TRUE),
control.compute = list(dic = TRUE, waic = TRUE, cpo = TRUE, config = TRUE)
)
# Smoothed rates (per person-year)
rate_mean <- fit_rate$summary.fitted.values$mean # rate
rate_lcl <- fit_rate$summary.fitted.values$`0.025quant`
rate_ucl <- fit_rate$summary.fitted.values$`0.975quant`
# If you prefer per 1,000 PY:
per <- 100
rate_mean_100 <- rate_mean * per
rate_lcl_100 <- rate_lcl * per
rate_ucl_100 <- rate_ucl * per
# --- 5B) BYM for excess risk and excess deaths ---
# Y ~ Poisson(E × RR_i); fitted values are RR
fit_rr <- INLA::inla(
Y ~ 1 + f(region_id, model = "bym2", graph = adj_file, scale.model = TRUE, hyper = hyper),
data = bym_df,
family = "poisson",
E = E, # expected deaths from reference
control.predictor = list(compute = TRUE),
control.compute = list(dic = TRUE, waic = TRUE, cpo = TRUE, config = TRUE)
)
RR_mean <- fit_rr$summary.fitted.values$mean
RR_lcl <- fit_rr$summary.fitted.values$`0.025quant`
RR_ucl <- fit_rr$summary.fitted.values$`0.975quant`
# Smoothed excess deaths = (RR - 1) * E
Excess_mean <- (RR_mean - 1) * bym_df$E
Excess_lcl <- (RR_lcl - 1) * bym_df$E
Excess_ucl <- (RR_ucl - 1) * bym_df$E
# Exceedance probability: P(RR > 1) = P(eta > 0)
if (!is.null(fit_rr$marginals.linear.predictor) &&
length(fit_rr$marginals.linear.predictor) == nrow(bym_df)) {
p_exceed_1 <- sapply(fit_rr$marginals.linear.predictor, function(m) 1 - INLA::inla.pmarginal(0, m))
} else {
# Monte Carlo backup
smp <- INLA::inla.posterior.sample(4000, fit_rr)
pred_names <- grep("^Predictor", rownames(smp[[1]]$latent), value = TRUE)
lp_mat <- sapply(smp, function(s) as.numeric(s$latent[pred_names]))
p_exceed_1 <- rowMeans(lp_mat > 0, na.rm = TRUE)
}
ps_pois <- inla.posterior.sample(2000, fit_rr, selection = list(Predictor = 1:nrow(fit_rr$summary.linear.predictor)))
eta_mat_pois <- sapply(ps_pois, function(s) s$latent[grep("^Predictor", rownames(s$latent)), , drop = TRUE])
p_exceed_1_mc <- rowMeans(eta_mat_pois > 0)
# --- 6) Bind results and map ---
res_cc1 <- bym_df|>
mutate(
rate_mean = rate_mean,
rate_lcl = rate_lcl,
rate_ucl = rate_ucl,
rate_mean_100 = rate_mean_100,
rate_lcl_100 = rate_lcl_100,
rate_ucl_100 = rate_ucl_100,
RR_mean = RR_mean,
RR_lcl = RR_lcl,
RR_ucl = RR_ucl,
Excess_mean = Excess_mean,
Excess_lcl = Excess_lcl,
Excess_ucl = Excess_ucl,
p_exceed_1 = p_exceed_1_mc
)
result_sf <- chile_sf|>
left_join(res_cc1, by = c("CC_1"="CC1"))
orden_geografico_romano <- c(
"XV", # Arica y Parinacota
"I", # Tarapacá
"II", # Antofagasta
"III", # Atacama
"IV", # Coquimbo
"V", # Valparaíso
"RM", # Metropolitana
"VI", # O'Higgins
"VII", # Maule
"XVI", # Ñuble
"VIII", # Biobío
"IX", # La Araucanía
"XIV", # Los Ríos
"X", # Los Lagos
"XI", # Aysén (lo incluyo por si aparece en tus datos)
"XII" # Magallanes
)
res_cc1 |>
dplyr::select(-dplyr::any_of(c("rate_mean","rate_lcl","rate_ucl", "p_exceed_1")))|>
dplyr::mutate(CC1 = factor(CC1, levels = orden_geografico_romano), )|>
dplyr::arrange(CC1)|>
knitr::kable("markdown", digits=1, caption= "Resultados del modelo BYM2 por Región")| CC1 | E | PY | Y | region_id | rate_mean_100 | rate_lcl_100 | rate_ucl_100 | RR_mean | RR_lcl | RR_ucl | Excess_mean | Excess_lcl | Excess_ucl |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| XV | 1.8 | 427.2 | 7 | 15 | 2.1 | 1.2 | 3.4 | 6.0 | 3.4 | 8.5 | 8.9 | 4.4 | 13.4 |
| I | 0.8 | 186.5 | 3 | 1 | 2.4 | 1.2 | 4.1 | 6.4 | 3.7 | 9.3 | 4.3 | 2.1 | 6.6 |
| II | 1.6 | 413.9 | 12 | 2 | 2.9 | 1.8 | 4.4 | 7.2 | 4.9 | 10.2 | 10.0 | 6.3 | 14.9 |
| III | 0.4 | 96.1 | 9 | 3 | 4.7 | 2.5 | 8.9 | 9.1 | 5.8 | 16.7 | 3.3 | 2.0 | 6.4 |
| IV | 0.8 | 200.5 | 8 | 4 | 3.4 | 2.1 | 5.6 | 7.7 | 5.2 | 11.8 | 5.5 | 3.4 | 8.8 |
| V | 3.2 | 765.6 | 17 | 7 | 2.5 | 1.7 | 3.5 | 6.4 | 4.4 | 8.5 | 17.3 | 10.8 | 23.8 |
| RM | 3.9 | 1251.5 | 32 | 6 | 2.6 | 1.9 | 3.5 | 7.6 | 5.8 | 10.1 | 25.6 | 18.5 | 35.2 |
| VI | 1.3 | 306.3 | 5 | 8 | 2.4 | 1.4 | 3.8 | 6.3 | 3.8 | 8.7 | 7.1 | 3.7 | 10.4 |
| VII | 1.0 | 220.3 | 6 | 9 | 3.0 | 1.7 | 4.9 | 6.9 | 4.4 | 9.9 | 5.9 | 3.4 | 8.9 |
| XVI | 0.0 | 0.0 | 0 | 16 | 3.7 | 1.7 | 7.4 | 7.3 | 4.2 | 11.8 | 0.0 | 0.0 | 0.0 |
| VIII | 3.3 | 551.2 | 19 | 10 | 3.7 | 2.5 | 5.3 | 6.7 | 4.6 | 8.8 | 18.8 | 12.0 | 25.8 |
| IX | 0.7 | 76.0 | 7 | 5 | 5.9 | 3.3 | 10.5 | 7.9 | 5.3 | 12.2 | 5.0 | 3.1 | 8.0 |
| XIV | 1.3 | 195.2 | 11 | 14 | 5.6 | 3.4 | 8.5 | 7.9 | 5.4 | 11.5 | 9.1 | 5.9 | 13.9 |
| X | 1.6 | 219.6 | 18 | 11 | 7.0 | 4.5 | 10.5 | 9.0 | 6.3 | 13.5 | 12.4 | 8.2 | 19.5 |
| XI | 0.0 | 0.0 | 0 | 12 | 6.3 | 2.6 | 12.6 | 8.1 | 4.6 | 13.9 | 0.0 | 0.0 | 0.0 |
| XII | 2.1 | 276.6 | 18 | 13 | 6.3 | 4.1 | 9.4 | 8.1 | 5.8 | 11.6 | 15.0 | 10.0 | 22.4 |
Visualización de resultados
Code
# Print plots
cowplot::plot_grid(
p_tratados +
coord_sf(expand = FALSE) +
theme(
plot.margin = margin(0, 0, 0, 0),
legend.spacing = unit(0.1, "lines"),
legend.key.size = unit(0.6, "lines")
), p_deaths +
coord_sf(expand = FALSE) +
theme(
plot.margin = margin(0, 0, 0, 0),
legend.spacing = unit(0.1, "lines"),
legend.key.size = unit(0.6, "lines")
),
ggplot(result_sf) +
geom_sf(aes(fill = RR_mean), color = "grey60", size = 0.2) +
scale_fill_gradientn(colors = rev(c(g052709, dark_teal, c5fa9b3, d_edee4, blanco)), na.value = "#DEDEE4", name = "RME Suavizado") +
#scale_fill_viridis_c(option = "C", name = "Smoothed RR") +
labs(title = NULL,#"BYM2-smoothed excess risk (SISTRAT vs reference)",
subtitle = NULL)+#"RR = observed / expected (age×sex region-specific rates)") +
theme_void(base_size = 14) +
#xlim(c(80,66.4))+
xlim(c(-75, -66.4))+
theme(legend.position = "left",
legend.title = element_text(size=12.9),
legend.text = element_text(size=11.9))+
ggrepel::geom_label_repel(
data = dplyr::mutate(result_sf,
# 1. Calculates a point guaranteed to be within the polygon boundaries.
point_on_surface = sf::st_point_on_surface(geometry),
# 2. Extracts the X (Longitude) and Y (Latitude) coordinates from that point.
lon = sf::st_coordinates(point_on_surface)[, "X"],
lat = sf::st_coordinates(point_on_surface)[, "Y"]), # Usamos los datos con las coordenadas calculadas
aes(x = lon, y = lat, label = CC_1), # Coordenadas y texto de la etiqueta
size = 4.6, # Tamaño del texto
min.segment.length = 0.5, # Longitud mínima de la línea de conexión
box.padding = 0.3, # Espacio alrededor del texto dentro de la caja
point.padding = 0.5, # Distancia de la caja a la geometría
# Fondo de la etiqueta
fill = alpha("white", 0.7), # Color de fondo semi-transparente
max.overlaps = Inf, # Increase repulsion effort to avoid all overlaps
seed = 2125, # For reproducible results
color = "black" # Color del texto
),
nrow = 1, # o nrow = 1 según tu disposición
align = "v", # o "h" para horizontal
rel_heights = c(1, 1),# ajusta proporciones si necesario
axis = "tblr", # alinea ejes si aplica
labels = NULL, # evita espacio extra por etiquetas
label_size = 0, # si usas labels, reduce tamaño
label_x = 0, label_y = 1,# controla posición de etiquetas
vjust = .5, hjust = 0, # controla alineación vertical/horizontal
scale = c(0.99, .99, 1.03) # reduce tamaño relativo de cada plot
)Warning: There was 1 warning in stopifnot(). ℹ In argument: point_on_surface = sf::st_point_on_surface(geometry). Caused by warning in st_point_on_surface.sfc(): ! st_point_on_surface may not give correct results for longitude/latitude data
Code
if(grepl("ubuntu",pak::system_r_platform())){
ggsave(paste0("map_con_RR.png"), dpi = 600, scale=1, width = 12, height = 7.38)
} else{
ggsave(paste0("map_con_RR.png"), dpi = 600, scale=1, width = 12, height = 7.38)
}Resumen
Code
# 1) Orden final (Totales arriba). Elimina la línea con `sections_order <- rev(sections_order)`
sections_order <- c(
"Totales",
"Por causas (aj. por sexo y edad)",
"Edad (aj. por sexo)",
"Sexo (aj. por grupo de edad)"
)
capitalize_first_letter_robust <- function(x) {
# 1. Limpiar espacios al inicio y final
x_clean <- trimws(x)
# 2. Capitalizar la primera letra
first_letter <- toupper(substring(x_clean, 1, 1))
# 3. Tomar el resto de la cadena a partir del segundo caracter
rest_of_string <- substring(x_clean, 2)
return(paste0(first_letter, rest_of_string))
}
# --- (después de tu asignación manual de `smr_resumen$section <- c(...)`) ---
# 3) Mantener solo esas secciones y fijar el orden como factor
keep_sections <- sections_order
smr_resumen2 <- smr_resumen |>
dplyr::filter(section %in% keep_sections) |>
dplyr::mutate(section = factor(section, levels = sections_order))|>
dplyr::filter(!res %in% c("Total ajustado por grupos de edad, sexo y año calendario", "Total ajustado por edad, sexo y año calendario"))|>
dplyr::mutate(section= dplyr::case_when(grepl("causas",section)~"Por causas (aj. por sexo\ny grupos de edad)",
TRUE~section))|>
dplyr::mutate(section= gsub("grupo de","grupos de", section))|>
dplyr::mutate(label = gsub("^Total", "", label))|>
dplyr::mutate(label = sapply(label, capitalize_first_letter_robust))
sections_order <- c(
"Totales",
"Por causas (aj. por sexo\ny grupos de edad)",
"Edad (aj. por sexo)",
"Sexo (aj. por grupos de edad)"
)
# Create display
display_list <- list()
y <- 0
for (sec in sections_order) {
header <- tibble(section = sec, label = NA_character_, obs = NA_real_, py = NA_real_, exp = NA_real_,
rme = NA_real_, lo = NA_real_, hi = NA_real_, is_header = TRUE, y = y)
display_list <- append(display_list, list(header))
y <- y + 1
items <- smr_resumen2 %>% filter(section == sec)
if (nrow(items) > 0) {
items <- items %>% mutate(is_header = FALSE, y = seq(y, by = 1, length.out = n()))
display_list <- append(display_list, list(items))
y <- y + nrow(items)
}
y <- y + 0.2
}
display <- bind_rows(display_list)
df_plot <- display %>% dplyr::filter(!is_header)
headers <- display %>% dplyr::filter(is_header) %>% dplyr::arrange(y)
# Colors
# 2) Paleta sincronizada con los nombres EXACTOS
section_cols <- c(
"Totales" = "#A17C6C80",
"Por causas (aj. por sexo y edad)" = "#906F0080",
"Edad (aj. por sexo)" = "#5FA9B380",
"Sexo (aj. por grupo de edad)" = "#D2A899"
)
# Calculate limits, start from 0.5
min_lo <- min(df_plot$lo, na.rm = TRUE)
max_hi <- max(df_plot$hi, na.rm = TRUE)
xmin <- 0.5 # Force start at 0.5
xmax <- exp(log(max_hi) + 0.0 * (log(max_hi) - log(min_lo)))
# Rects
rects <- headers %>%
transmute(
section,
ymin = y + 0.4,
ymax = lead(y) - 0.4
)
rects$ymax[is.na(rects$ymax)] <- max(display$y, na.rm = TRUE) + 0.6
rects <- rects %>%
mutate(xmin = xmin * 1.005,
xmax = xmax / 1.005)
# Labels position (moved left)
x_left_lab <- xmin * 0.5 # Moved further left
cap_h <- 0.25
p2 <- ggplot() +
geom_rect(data = rects,
aes(xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax),
inherit.aes = FALSE, fill = "#f7f7f7", color = NA, alpha = 0.6) +
geom_vline(xintercept = 1, linetype = "dashed", color = "#444444", linewidth = 1) +
geom_segment(data = df_plot,
aes(x = lo, xend = hi, y = y, yend = y, color = section),
linewidth = 0.9) +
geom_segment(data = df_plot,
aes(x = lo, xend = lo, y = y - cap_h, yend = y + cap_h, color = section),
linewidth = 0.9) +
geom_segment(data = df_plot,
aes(x = hi, xend = hi, y = y - cap_h, yend = y + cap_h, color = section),
linewidth = 0.9) +
geom_point(data = df_plot,
aes(x = rme, y = y, fill = section, color = section),
shape = 22, size=3,stroke = 0.6) +
# Label over point with RME to 1 decimal, larger size
geom_text(data = df_plot,
aes(x = rme, y = y + 0.3, label = sprintf("%.1f", rme)),
size = 5, color = "black", hjust = 0.5, vjust = 0) +
geom_text(data = df_plot,
aes(x = x_left_lab-.25, y = y, label = paste0(" ", label)),
hjust = 0, vjust = 0.5, size = 3.6, color = "#222222", inherit.aes = FALSE) +
geom_text(data = left_join(headers, rects[, c("section", "ymax")], by = "section") %>% mutate(y = ymax + .4),
aes(x = x_left_lab-.25, y = y, label = section),
hjust = 0, vjust = 0.5, fontface = "bold", size = 4.0, inherit.aes = FALSE) +
scale_color_manual(values = section_cols, name = "Sección") +
scale_fill_manual(values = section_cols, guide = "none") +
#scale_size(range = c(2.6, 6.2), guide = "none") +
scale_x_continuous(
trans = scales::log_trans(base = exp(1)),
breaks = c(0.5, 1, 2, 4, 8, 16),
labels = function(x) format(x, big.mark = ".", decimal.mark = ",", trim = TRUE),
expand = ggplot2::expansion(mult = c(0.3, 0.0)) # No expansion
) +
coord_cartesian(xlim = c(0.5, xmax)) + # Force xlim from 0.5
labs(
title = NULL,
subtitle = NULL,
x = "RME (IC 95%), escala logarítmica",
y = NULL
) +
theme_minimal(base_size = 15) +
theme(
panel.grid.major.y = element_blank(),
panel.grid.minor = element_blank(),
axis.text.y = element_blank(),
axis.title.y = element_blank(),
plot.title.position = "plot",
legend.position = "none"
)
print(p2)Warning in scale_x_continuous(trans = scales::log_trans(base = exp(1)), : log-2.718282 transformation introduced infinite values. log-2.718282 transformation introduced infinite values.
Code
if(grepl("ubuntu",pak::system_r_platform())){
ggsave(paste0("forestplot_RME_log_corr.png"), dpi = 600)
} else{
ggsave(paste0("forestplot_RME_log_corr.png"), dpi = 600)
}Warning in scale_x_continuous(trans = scales::log_trans(base = exp(1)), : log-2.718282 transformation introduced infinite values. log-2.718282 transformation introduced infinite values.
Session info
Code
#|echo: true
#|error: true
#|message: true
#|paged.print: true
message(paste0("R library: ", Sys.getenv("R_LIBS_USER")))Code
message(paste0("Date: ",withr::with_locale(new = c('LC_TIME' = 'C'), code =Sys.time())))Code
message(paste0("Editor context: ", path))Error: objeto ‘path’ no encontrado
Code
cat("quarto version: "); quarto::quarto_version()quarto version: [1] '1.7.32'Code
sesion_info <- devtools::session_info()
dplyr::select(
tibble::as_tibble(sesion_info$packages),
c(package, loadedversion, source)
)|>
DT::datatable(filter = 'top', colnames = c('Row number' =1,'Package' = 2, 'Version'= 3),
caption = htmltools::tags$caption(
style = 'caption-side: top; text-align: left;',
'', htmltools::em('R packages')),
options=list(
initComplete = htmlwidgets::JS(
"function(settings, json) {",
"$(this.api().tables().body()).css({
'font-family': 'Helvetica Neue',
'font-size': '70%',
'code-inline-font-size': '15%',
'white-space': 'nowrap',
'line-height': '0.75em',
'min-height': '0.5em'
});",
"}")))Code
#|echo: true
#|error: true
#|message: true
#|paged.print: true
#|class-output: center-table
reticulate::py_list_packages()|>
DT::datatable(filter = 'top', colnames = c('Row number' =1,'Package' = 2, 'Version'= 3),
caption = htmltools::tags$caption(
style = 'caption-side: top; text-align: left;',
'', htmltools::em('Python packages')),
options=list(
initComplete = htmlwidgets::JS(
"function(settings, json) {",
"$(this.api().tables().body()).css({
'font-family': 'Helvetica Neue',
'font-size': '70%',
'code-inline-font-size': '15%',
'white-space': 'nowrap',
'line-height': '0.75em',
'min-height': '0.5em'
});",
"}"))) Save
Code
if(grepl("ubuntu",pak::system_r_platform())){
if (rstudioapi::isAvailable()) {
# Code that needs RStudio
wdpath<- dirname(dirname(rstudioapi::documentPath()))
} else {
# Code for non-RStudio environments (e.g., command line, server)
# This part should use a portable method like knitr::current_input() or getwd()
wdpath <- dirname(resolve_doc_dir())
}
}
save.image(paste0(wdpath,"/congresosp2.RData"))
cat(paste0("Espacio de trabajo guardado en ",wdpath,"/congresosp2.RData\n"))
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
cat("Copy renv lock into cons folder\n")
if (Sys.getenv("RSTUDIO_SESSION_TYPE") == "server" || file.exists("/.dockerenv")) {
message("Running on RStudio Server or inside Docker. Folder copy skipped.")
} else {
source_folder <-
destination_folder <- paste0(wdpath,"cons/renv")
# Copy the folder recursively
file.copy(paste0(wdpath,"renv.lock"), paste0(wdpath,"cons/renv.lock"), overwrite = TRUE)
message("Renv lock copy performed.")
}Code
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
time_after_dedup2<-Sys.time()
paste0("Time in markdown: ");time_after_dedup2-time_before_dedup2Espacio de trabajo guardado en /home/andres/Escritorio/SER2025/congresosp2.RData
Copy renv lock into cons folder
[1] "Time in markdown: "
Time difference of 16.21523 days